<!DOCTYPE html>
<html>
<head>
  <meta charset="utf-8">
  
  <title>一个具体能量方程的解析 | Giskard&#39;s CFD Learning Tricks</title>
  <meta name="viewport" content="width=device-width, initial-scale=1, maximum-scale=1">
  <meta name="description" content="本篇来看一个具体的能量方程，以 twoPhaseEulerFoam 的 EEqn.H 为例。">
<meta property="og:type" content="article">
<meta property="og:title" content="一个具体能量方程的解析">
<meta property="og:url" content="http://xiaopingqiu.github.io/2016/06/25/thermophysics4/index.html">
<meta property="og:site_name" content="Giskard's CFD Learning Tricks">
<meta property="og:description" content="本篇来看一个具体的能量方程，以 twoPhaseEulerFoam 的 EEqn.H 为例。">
<meta name="twitter:card" content="summary">
<meta name="twitter:title" content="一个具体能量方程的解析">
<meta name="twitter:description" content="本篇来看一个具体的能量方程，以 twoPhaseEulerFoam 的 EEqn.H 为例。">
  
  
    <link rel="icon" href="/favicon.png">
  
  <link href="//fonts.googleapis.com/css?family=Source+Code+Pro" rel="stylesheet" type="text/css">
  <link rel="stylesheet" href="/css/style.css" type="text/css">
  
<!-- Google Analytics -->
<script type="text/javascript">
(function(i,s,o,g,r,a,m){i['GoogleAnalyticsObject']=r;i[r]=i[r]||function(){
(i[r].q=i[r].q||[]).push(arguments)},i[r].l=1*new Date();a=s.createElement(o),
m=s.getElementsByTagName(o)[0];a.async=1;a.src=g;m.parentNode.insertBefore(a,m)
})(window,document,'script','//www.google-analytics.com/analytics.js','ga');

##ga('create', '[object Object]', 'auto');
ga('create', 'UA-62501686-1', 'auto');
ga('send', 'pageview');

</script>
<!-- End Google Analytics -->


</head>
<body>
  <div id="container">
    <div id="wrap">
      <header id="header">
  <div id="banner"></div>
  <div id="header-outer" class="outer">
    <div id="header-title" class="inner">
      <h1 id="logo-wrap">
        <a href="/" id="logo">Giskard&#39;s CFD Learning Tricks</a>
      </h1>
      
        <h2 id="subtitle-wrap">
          <a href="/" id="subtitle">CFD and Scientific Computing</a>
        </h2>
      
    </div>
    <div id="header-inner" class="inner">
      <nav id="main-nav">
        <a id="main-nav-toggle" class="nav-icon"></a>
        
          <a class="main-nav-link" href="/">Home</a>
        
          <a class="main-nav-link" href="/archives">Archives</a>
        
          <a class="main-nav-link" href="/atom.xml">Rss</a>
        
          <a class="main-nav-link" href="/about">About</a>
        
      </nav>
      <nav id="sub-nav">
        
        <a id="nav-search-btn" class="nav-icon" title="Search"></a>
      </nav>
      <div id="search-form-wrap">
        <form action="//google.com/search" method="get" accept-charset="UTF-8" class="search-form"><input type="search" name="q" results="0" class="search-form-input" placeholder="Search"><button type="submit" class="search-form-submit">&#xF002;</button><input type="hidden" name="q" value="site:http://xiaopingqiu.github.io"></form>
      </div>
    </div>
  </div>
</header>

      <div class="outer">
        <section id="main"><article id="post-thermophysics4" class="article article-type-post" itemscope itemprop="blogPost">
  <div class="article-meta">
    <a href="/2016/06/25/thermophysics4/" class="article-date">
  <time datetime="2016-06-25T07:48:02.000Z" itemprop="datePublished">2016-06-25</time>
</a>
    
  <div class="article-category">
    <a class="article-category-link" href="/categories/OpenFOAM/">OpenFOAM</a>
  </div>

  </div>
  <div class="article-inner">
    
    
      <header class="article-header">
        
  
    <h1 class="article-title" itemprop="name">
      一个具体能量方程的解析
    </h1>
  

      </header>
    
    <div class="article-entry" itemprop="articleBody">
      
        <p>本篇来看一个具体的能量方程，以 <code>twoPhaseEulerFoam</code> 的 <code>EEqn.H</code> 为例。</p>
<a id="more"></a>
<figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br><span class="line">26</span><br><span class="line">27</span><br><span class="line">28</span><br><span class="line">29</span><br><span class="line">30</span><br><span class="line">31</span><br><span class="line">32</span><br><span class="line">33</span><br><span class="line">34</span><br><span class="line">35</span><br><span class="line">36</span><br><span class="line">37</span><br><span class="line">38</span><br><span class="line">39</span><br></pre></td><td class="code"><pre><span class="line">&#123;</span><br><span class="line">    volScalarField&amp; he1 = thermo1.he();</span><br><span class="line">    volScalarField&amp; he2 = thermo2.he();</span><br><span class="line"></span><br><span class="line">    <span class="function">volScalarField <span class="title">Cpv1</span><span class="params">("Cpv1", thermo1.Cpv()</span>)</span>;</span><br><span class="line">    <span class="function">volScalarField <span class="title">Cpv2</span><span class="params">("Cpv2", thermo2.Cpv()</span>)</span>;</span><br><span class="line"></span><br><span class="line">    <span class="function">volScalarField <span class="title">heatTransferCoeff</span><span class="params">(fluid.heatTransferCoeff()</span>)</span>;</span><br><span class="line"></span><br><span class="line">    <span class="function">fvScalarMatrix <span class="title">he1Eqn</span></span><br><span class="line">    <span class="params">(</span><br><span class="line">        fvm::ddt(alpha1, rho1, he1)</span> + fvm::<span class="title">div</span><span class="params">(alphaRhoPhi1, he1)</span></span><br><span class="line">      - fvm::<span class="title">Sp</span><span class="params">(contErr1, he1)</span></span><br><span class="line"></span><br><span class="line">      + fvc::<span class="title">ddt</span><span class="params">(alpha1, rho1, K1)</span> + fvc::<span class="title">div</span><span class="params">(alphaRhoPhi1, K1)</span></span><br><span class="line">      - contErr1*K1</span><br><span class="line">      + <span class="params">(</span><br><span class="line">            he1.name()</span> </span>== thermo1.phasePropertyName(<span class="string">"e"</span>)</span><br><span class="line">          ? fvc::ddt(alpha1)*p + fvc::div(alphaPhi1, p)</span><br><span class="line">          : -alpha1*dpdt</span><br><span class="line">        )</span><br><span class="line"></span><br><span class="line">      - fvm::laplacian</span><br><span class="line">        (</span><br><span class="line">            fvc::interpolate(alpha1)</span><br><span class="line">           *fvc::interpolate(thermo1.alphaEff(phase1.turbulence().mut())),</span><br><span class="line">            he1</span><br><span class="line">        )</span><br><span class="line">    );</span><br><span class="line"></span><br><span class="line">    he1Eqn.relax();</span><br><span class="line"></span><br><span class="line">    he1Eqn -=</span><br><span class="line">    (</span><br><span class="line">        heatTransferCoeff*(thermo2.T() - thermo1.T())</span><br><span class="line">      + heatTransferCoeff*he1/Cpv1</span><br><span class="line">      - fvm::Sp(heatTransferCoeff/Cpv1, he1)</span><br><span class="line">      + fvOptions(alpha1, rho1, he1)</span><br><span class="line">    );</span><br></pre></td></tr></table></figure>
<p>对应的能量方程为（忽略fvOptions）<br>$$<br>\alpha \rho \frac{\partial (\mathrm{he})}{\partial t} + \alpha \rho U\cdot \nabla(\mathrm{he}) + \alpha \rho \frac{\partial (\mathrm{K})}{\partial t} + \alpha \rho U\cdot \nabla\mathrm{K} + \\<br>\begin{cases}<br>p\cdot\dfrac{\partial \alpha}{\partial t} + \nabla \cdot (\alpha U p) , &amp; \mbox{ if } he.name == \mbox{“e”} \\<br>-\alpha \dfrac{\partial p}{\partial t}, &amp; \mbox{  if } he.name == \mbox{“h”}<br>\end{cases} \\<br>-\nabla \cdot \big(\alpha \cdot \alpha_{eff} \nabla (\mathrm{he}) \big) - \gamma(T_2 - T_1) = 0<br>$$</p>
<p>代码里剩下的两项，<br><figure class="highlight stylus"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br></pre></td><td class="code"><pre><span class="line">+ heatTransferCoeff*he1/Cpv1</span><br><span class="line">- fvm::<span class="function"><span class="title">Sp</span><span class="params">(heatTransferCoeff/Cpv1, he1)</span></span></span><br></pre></td></tr></table></figure></p>
<p>含义暂不明。这两项，其实是同一个公式，只是前者是显示处理，后者用了隐式源项，估计是为了数值稳定性的目的而构建的。</p>
<p>前面提过，对于如下设置，<br><figure class="highlight nginx"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br></pre></td><td class="code"><pre><span class="line"><span class="title">thermoType</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="title">type</span>             heRhoThermo;</span><br><span class="line">    <span class="title">mixture</span>          pureMixture;</span><br><span class="line">    <span class="title">transport</span>        const;</span><br><span class="line">    <span class="title">thermo</span>           hConst;</span><br><span class="line">    <span class="title">equationOfState</span>  perfectGas;</span><br><span class="line">    <span class="title">specie</span>           specie;</span><br><span class="line">    <span class="title">energy</span>           sensibleInternalEnergy;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>最终，<code>thermo</code> 指针指向的是 <code>heRhoThermo</code> 类的对象，所以，从 <code>heRhoThermo</code> 类的构造函数看起：<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">BasicPsiThermo</span>, class <span class="type">MixtureType</span>&gt;</span><br><span class="line"><span class="type">Foam</span>::heRhoThermo&lt;<span class="type">BasicPsiThermo</span>, <span class="type">MixtureType</span>&gt;::heRhoThermo</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> fvMesh&amp; mesh,</span><br><span class="line">    <span class="keyword">const</span> word&amp; phaseName</span><br><span class="line">)</span><br><span class="line">:</span><br><span class="line">    heThermo&lt;<span class="type">BasicPsiThermo</span>, <span class="type">MixtureType</span>&gt;(mesh, phaseName)</span><br><span class="line">&#123;</span><br><span class="line">    calculate(); // 构造函数调用 calculate 函数来初始化所有的热物理相关量</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>可见，构造函数里调用了 <code>calculate</code> 函数，前面提过，这个函数的作用是更新各个热物理相关量。</p>
<p>接下来一个一个来看里面涉及到的函数。</p>
<h5 id="he"><code>he</code></h5><p><code>he</code> 其实是 “h or e”，具体是焓，还是内能，取决于 <code>energy</code> 这一项的设置。 <code>he</code> 函数在 <code>heThermo</code> 类中定义，返回的是数据成员 <code>he_</code>，所以这里需要看一下数据成员 <code>he_</code> 的初始化：<br><figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br><span class="line">26</span><br><span class="line">27</span><br><span class="line">28</span><br><span class="line">29</span><br><span class="line">30</span><br><span class="line">31</span><br><span class="line">32</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;<span class="keyword">class</span> BasicThermo, <span class="keyword">class</span> MixtureType&gt;</span><br><span class="line">Foam::heThermo&lt;BasicThermo, MixtureType&gt;::heThermo</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> fvMesh&amp; mesh,</span><br><span class="line">    <span class="keyword">const</span> dictionary&amp; dict,</span><br><span class="line">    <span class="keyword">const</span> word&amp; phaseName</span><br><span class="line">)</span><br><span class="line">:</span><br><span class="line">    BasicThermo(mesh, dict, phaseName),</span><br><span class="line">    MixtureType(*<span class="keyword">this</span>, mesh),</span><br><span class="line"></span><br><span class="line">    he_</span><br><span class="line">    (</span><br><span class="line">        IOobject</span><br><span class="line">        (</span><br><span class="line">            BasicThermo::phasePropertyName</span><br><span class="line">            (</span><br><span class="line">                MixtureType::thermoType::heName()</span><br><span class="line">            ),</span><br><span class="line">            mesh.time().timeName(),</span><br><span class="line">            mesh,</span><br><span class="line">            IOobject::NO_READ,</span><br><span class="line">            IOobject::NO_WRITE</span><br><span class="line">        ),</span><br><span class="line">        mesh,</span><br><span class="line">        dimEnergy/dimMass,</span><br><span class="line">        <span class="keyword">this</span>-&gt;heBoundaryTypes(),</span><br><span class="line">        <span class="keyword">this</span>-&gt;heBoundaryBaseTypes()</span><br><span class="line">    )</span><br><span class="line">&#123;</span><br><span class="line">    init();</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这里调用的 <code>init</code> 函数的内容为<br><figure class="highlight ocaml"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br></pre></td><td class="code"><pre><span class="line">template&lt;<span class="keyword">class</span> <span class="type">BasicThermo</span>, <span class="keyword">class</span> <span class="type">MixtureType</span>&gt;</span><br><span class="line">void <span class="type">Foam</span>::heThermo&lt;<span class="type">BasicThermo</span>, <span class="type">MixtureType</span>&gt;::init<span class="literal">()</span></span><br><span class="line">&#123;</span><br><span class="line">    scalarField&amp; heCells = he_.internalField<span class="literal">()</span>;</span><br><span class="line">    const scalarField&amp; pCells = this-&gt;p_.internalField<span class="literal">()</span>;</span><br><span class="line">    const scalarField&amp; <span class="type">TCells</span> = this-&gt;<span class="type">T_</span>.internalField<span class="literal">()</span>;</span><br><span class="line"></span><br><span class="line">    forAll(heCells, celli)</span><br><span class="line">    &#123;</span><br><span class="line">        heCells[celli] =</span><br><span class="line">            this-&gt;cellMixture(celli).<span class="type">HE</span>(pCells[celli], <span class="type">TCells</span>[celli]);</span><br><span class="line">    &#125;</span><br><span class="line"></span><br><span class="line">    forAll(he_.boundaryField<span class="literal">()</span>, patchi)</span><br><span class="line">    &#123;</span><br><span class="line">        he_.boundaryField<span class="literal">()</span>[patchi] == he</span><br><span class="line">        (</span><br><span class="line">            this-&gt;p_.boundaryField<span class="literal">()</span>[patchi],</span><br><span class="line">            this-&gt;<span class="type">T_</span>.boundaryField<span class="literal">()</span>[patchi],</span><br><span class="line">            patchi</span><br><span class="line">        );</span><br><span class="line">    &#125;</span><br><span class="line"></span><br><span class="line">    this-&gt;heBoundaryCorrection(he_);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这里调用了 <code>HE</code> 函数来初始化 <code>he_</code> 的内部场，并对调用另一个三个数的 <code>he</code> 函数其边界条件进行了修正：<br><figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;<span class="keyword">class</span> BasicThermo, <span class="keyword">class</span> MixtureType&gt;</span><br><span class="line">Foam::tmp&lt;Foam::scalarField&gt; Foam::heThermo&lt;BasicThermo, MixtureType&gt;::he</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> scalarField&amp; p,</span><br><span class="line">    <span class="keyword">const</span> scalarField&amp; T,</span><br><span class="line">    <span class="keyword">const</span> label patchi</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    tmp&lt;scalarField&gt; the(<span class="keyword">new</span> scalarField(T.size()));</span><br><span class="line">    scalarField&amp; he = the();</span><br><span class="line"></span><br><span class="line">    forAll(T, facei)</span><br><span class="line">    &#123;</span><br><span class="line">        he[facei] =</span><br><span class="line">            <span class="keyword">this</span>-&gt;patchFaceMixture(patchi, facei).HE(p[facei], T[facei]);</span><br><span class="line">            <span class="comment">// 本质上还是调用了 HE 函数</span></span><br><span class="line">    &#125;</span><br><span class="line">    <span class="keyword">return</span> the;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>再来看 <code>HE</code> 函数，这个函数看名字和参数，应该是根据压力和温度来计算能量的，其定义在 <code>species::thermo</code> 类：<br><figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;<span class="keyword">class</span> Thermo, <span class="keyword">template</span>&lt;<span class="keyword">class</span>&gt; <span class="keyword">class</span> Type&gt;</span><br><span class="line"><span class="keyword">inline</span> Foam::scalar</span><br><span class="line">Foam::species::thermo&lt;Thermo, Type&gt;::HE(<span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> Type&lt;thermo&lt;Thermo, Type&gt; &gt;::HE(*<span class="keyword">this</span>, p, T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这里，由于能量最终是什么形式，取决于 <code>energy</code> 关键字对应的类，所以，这里也是调用了定义在前面提到的 <code>energy variable</code> 类中的 <code>HE</code> 函数，以 <code>sensibleInternalEnergy</code> 为例：<br><figure class="highlight stata"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">scalar</span> <span class="keyword">HE</span></span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> Thermo&amp; thermo,</span><br><span class="line">    <span class="keyword">const</span> <span class="keyword">scalar</span> p,</span><br><span class="line">    <span class="keyword">const</span> <span class="keyword">scalar</span> T</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> thermo.Es(p, T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>可见，其返回的是 <code>species::thermo</code> 类的 <code>Es</code> 函数，<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar</span><br><span class="line"><span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::<span class="type">Es</span>(<span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> this-&gt;es(p, T)/this-&gt;W();</span><br><span class="line">&#125;</span><br><span class="line"></span><br><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar</span><br><span class="line"><span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::es(<span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> this-&gt;hs(p, T) - p*this-&gt;W()/this-&gt;rho(p, T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p> <code>hs</code> 函数定义在 <code>thermo</code> 类型的类中，以 <code>hConstThermo</code> 类为例：<br><figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;<span class="keyword">class</span> equationOfState&gt;</span><br><span class="line"><span class="keyword">inline</span> Foam::scalar Foam::hConstThermo&lt;equationOfState&gt;::hs</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> Cp_*T;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p><code>hs</code> 表示的是显焓，等于 <code>Cp_*T</code> 。 <code>es</code> 是内能，根据焓的定义，$H=U+pV$。代码中的 <code>hs</code> 和 <code>es</code> 都是 <code>J/kMol</code> 的量纲，所以，$es=hs-pV/n$ 。以理想气体状态方程为例，$pV=nRT$，或者写成 $pM=\rho RT$，得 $pV/n = RT = pM/\rho$ 。</p>
<p>注意，这里的 <code>Cp_</code>，在字典文件里给的是 <code>J/(kg.K)</code> 量纲的，但是在构造函数中，将其转成了 <code>J/(kmol.K)</code> 的量纲：<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class equationOfState&gt;</span><br><span class="line"><span class="type">Foam</span>::hConstThermo&lt;equationOfState&gt;::hConstThermo(<span class="keyword">const</span> dictionary&amp; dict)</span><br><span class="line">:</span><br><span class="line">    equationOfState(dict),</span><br><span class="line">    <span class="type">Cp_</span>(readScalar(dict.subDict(<span class="string">"thermodynamics"</span>).lookup(<span class="string">"Cp"</span>))),</span><br><span class="line">    <span class="type">Hf_</span>(readScalar(dict.subDict(<span class="string">"thermodynamics"</span>).lookup(<span class="string">"Hf"</span>)))</span><br><span class="line">&#123;</span><br><span class="line">    <span class="type">Cp_</span> *= this-&gt;W();</span><br><span class="line">    <span class="type">Hf_</span> *= this-&gt;W();</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>所以，<code>hs</code>， <code>es</code> 是 <code>J/kmol</code> ； <code>Es</code>， <code>HE</code> 是 <code>J/kg</code>。</p>
<h5 id="Cpv">Cpv</h5><p>这个函数定义在 <code>heThermo</code> 类中。<br><figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br><span class="line">26</span><br><span class="line">27</span><br><span class="line">28</span><br><span class="line">29</span><br><span class="line">30</span><br><span class="line">31</span><br><span class="line">32</span><br><span class="line">33</span><br><span class="line">34</span><br><span class="line">35</span><br><span class="line">36</span><br><span class="line">37</span><br><span class="line">38</span><br><span class="line">39</span><br><span class="line">40</span><br><span class="line">41</span><br><span class="line">42</span><br><span class="line">43</span><br><span class="line">44</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;<span class="keyword">class</span> BasicThermo, <span class="keyword">class</span> MixtureType&gt;</span><br><span class="line">Foam::tmp&lt;Foam::volScalarField&gt;</span><br><span class="line">Foam::heThermo&lt;BasicThermo, MixtureType&gt;::Cpv() <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">const</span> fvMesh&amp; mesh = <span class="keyword">this</span>-&gt;T_.mesh();</span><br><span class="line"></span><br><span class="line">    tmp&lt;volScalarField&gt; tCpv</span><br><span class="line">    (</span><br><span class="line">        <span class="keyword">new</span> volScalarField</span><br><span class="line">        (</span><br><span class="line">            IOobject</span><br><span class="line">            (</span><br><span class="line">                <span class="string">"Cpv"</span>,</span><br><span class="line">                mesh.time().timeName(),</span><br><span class="line">                mesh,</span><br><span class="line">                IOobject::NO_READ,</span><br><span class="line">                IOobject::NO_WRITE</span><br><span class="line">            ),</span><br><span class="line">            mesh,</span><br><span class="line">            dimEnergy/dimMass/dimTemperature</span><br><span class="line">        )</span><br><span class="line">    );</span><br><span class="line">    volScalarField&amp; cpv = tCpv();</span><br><span class="line"></span><br><span class="line">    forAll(<span class="keyword">this</span>-&gt;T_, celli)</span><br><span class="line">    &#123;</span><br><span class="line">        cpv[celli] =</span><br><span class="line">            <span class="keyword">this</span>-&gt;cellMixture(celli).Cpv(<span class="keyword">this</span>-&gt;p_[celli], <span class="keyword">this</span>-&gt;T_[celli]);</span><br><span class="line">    &#125;</span><br><span class="line">    forAll(<span class="keyword">this</span>-&gt;T_.boundaryField(), patchi)</span><br><span class="line">    &#123;</span><br><span class="line">        <span class="keyword">const</span> fvPatchScalarField&amp; pp = <span class="keyword">this</span>-&gt;p_.boundaryField()[patchi];</span><br><span class="line">        <span class="keyword">const</span> fvPatchScalarField&amp; pT = <span class="keyword">this</span>-&gt;T_.boundaryField()[patchi];</span><br><span class="line">        fvPatchScalarField&amp; pCpv = cpv.boundaryField()[patchi];</span><br><span class="line"></span><br><span class="line">        forAll(pT, facei)</span><br><span class="line">        &#123;</span><br><span class="line">            pCpv[facei] =</span><br><span class="line">                <span class="keyword">this</span>-&gt;patchFaceMixture(patchi, facei).Cpv(pp[facei], pT[facei]);</span><br><span class="line">        &#125;</span><br><span class="line">    &#125;</span><br><span class="line"></span><br><span class="line">    <span class="keyword">return</span> tCpv;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这个函数，创建了一个 <code>tmp&lt;volScalarField&gt;</code>，然后调用定义在 <code>species::thermo</code> 类中的两参数 <code>Cpv</code> 函数来对场量进行初始化，这个函数的形式如下<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar</span><br><span class="line"><span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::<span class="type">Cpv</span>(<span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> this-&gt;cpv(p, T)/this-&gt;W();</span><br><span class="line">&#125;</span><br><span class="line"></span><br><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar</span><br><span class="line"><span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::cpv(<span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> <span class="type">Type</span>&lt;thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt; &gt;::cpv(*this, p, T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p><code>cpv</code> 函数返回的是定义在 <code>energy variable</code> 类中的三参数 <code>cpv</code> 函数，对于 <code>sensibleInternalEnergy</code>，<br><figure class="highlight stata"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">scalar</span> cpv</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> Thermo&amp; thermo,</span><br><span class="line">    <span class="keyword">const</span> <span class="keyword">scalar</span> p,</span><br><span class="line">    <span class="keyword">const</span> <span class="keyword">scalar</span> T</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> thermo.cv(p, T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这里返回的是 <code>species::thermo</code> 类的 <code>cv</code> 函数，<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar</span><br><span class="line"><span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::cv(<span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> this-&gt;cp(p, T) - this-&gt;cpMcv(p, T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这里的 <code>cp</code> 函数定义在 <code>thermo</code> 类型的类中，以 <code>hConst</code> 为例<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class equationOfState&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar <span class="type">Foam</span>::hConstThermo&lt;equationOfState&gt;::cp</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> scalar p,</span><br><span class="line">    <span class="keyword">const</span> scalar T</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> <span class="type">Cp_</span>; // 量纲是 J/(kmol.K)</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p><code>cpMcv</code> 的定义在状态方程类中，以 <code>perfectGas</code> 为例<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Specie</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar <span class="type">Foam</span>::perfectGas&lt;<span class="type">Specie</span>&gt;::cpMcv(scalar, scalar) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> this-&gt;<span class="type">RR</span>; // 量纲是 J/(kmol.K)，所以值应该是 <span class="number">8314</span></span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<h5 id="alphaEff">alphaEff</h5><p>这个函数需要一个参数，其定义在 <code>heThermo</code> 类中<br><figure class="highlight elixir"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br></pre></td><td class="code"><pre><span class="line">template&lt;class <span class="constant">BasicThermo,</span> class <span class="constant">MixtureType&gt;</span></span><br><span class="line"><span class="constant">Foam:</span><span class="symbol">:tmp&lt;Foam</span><span class="symbol">:</span><span class="symbol">:volScalarField&gt;</span></span><br><span class="line"><span class="constant">Foam:</span><span class="symbol">:heThermo&lt;BasicThermo</span>, <span class="constant">MixtureType&gt;</span><span class="symbol">:</span><span class="symbol">:alphaEff</span></span><br><span class="line">(</span><br><span class="line">    const volScalarField&amp; alphat</span><br><span class="line">) const</span><br><span class="line">&#123;</span><br><span class="line">    tmp&lt;<span class="constant">Foam:</span><span class="symbol">:volScalarField&gt;</span> alphaEff(this-&gt;<span class="constant">CpByCpv(</span>)*(this-&gt;alpha<span class="constant">_</span> + alphat));</span><br><span class="line">    alphaEff().rename(<span class="string">"alphaEff"</span>);</span><br><span class="line">    <span class="keyword">return</span> alphaEff;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这里， 无参数的 <code>CpByCpv</code> 函数定义在 <code>species::thermo</code> 类中，最终调用的是 <code>energy varialble</code> 类中的 <code>CpByCpv</code> 函数，如果是内能形式的，则返回 <code>thermo.gamma(p, T)</code> ，焓形式则返回 <code>1</code>。 <code>gamma</code> 的定义在 <code>species::thermo</code><br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar</span><br><span class="line"><span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::gamma(<span class="keyword">const</span> scalar p, <span class="keyword">const</span> scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    scalar cp = this-&gt;cp(p, T);</span><br><span class="line">    <span class="keyword">return</span> cp/(cp - this-&gt;cpMcv(p, T));</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p><code>alpha_</code> 是层流能量扩散系数，定义在 <code>basicThermo</code> 类，并在该类的构造函数中初始化为零。在 <code>heRhoThermo</code> 类的 <code>calculate</code> 函数中，对其进行了更新<br><figure class="highlight prolog"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br></pre></td><td class="code"><pre><span class="line"><span class="atom">scalarField</span>&amp; <span class="atom">alphaCells</span> = <span class="atom">this</span>-&gt;<span class="atom">alpha_</span>.<span class="atom">internalField</span>();</span><br><span class="line"><span class="atom">alphaCells</span>[<span class="atom">celli</span>] = <span class="atom">mixture_</span>.<span class="atom">alphah</span>(<span class="atom">pCells</span>[<span class="atom">celli</span>], <span class="name">TCells</span>[<span class="atom">celli</span>]);</span><br><span class="line"></span><br><span class="line"><span class="atom">fvPatchScalarField</span>&amp; <span class="atom">palpha</span> = <span class="atom">this</span>-&gt;<span class="atom">alpha_</span>.<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line"><span class="atom">palpha</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">alphah</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br></pre></td></tr></table></figure></p>
<p>可见， <code>alpha_</code> 的值是通过 <code>alphah</code> 函数来计算更新的，这个函数的定义在 <code>trasnport</code> 模型里，以 <code>constTransport</code> 为例<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar <span class="type">Foam</span>::constTransport&lt;<span class="type">Thermo</span>&gt;::alphah</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> scalar p,</span><br><span class="line">    <span class="keyword">const</span> scalar T</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> mu(p, T)*rPr_;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>返回粘度与普朗特数的比值。<br>至于 <code>alphat</code> ，则是函数 <code>alphaEff</code> 的参数，根据开头的代码可知， <code>alphat</code> 其实是 <code>mut</code> 。<br>只是，暂时不知道为什么有效热扩散系数 <code>alphaEff = CpByCpv * (alpha + alphat)</code>。</p>
<p>构建起能量方程后，就该对其进行求解了。<br><figure class="highlight ocaml"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br></pre></td><td class="code"><pre><span class="line">fvOptions.constrain(he1Eqn);</span><br><span class="line">   he1Eqn.solve<span class="literal">()</span>;</span><br><span class="line"></span><br><span class="line">   fvOptions.constrain(he2Eqn);</span><br><span class="line">   he2Eqn.solve<span class="literal">()</span>;</span><br><span class="line"></span><br><span class="line">   thermo1.correct<span class="literal">()</span>;</span><br><span class="line">   <span class="type">Info</span>&lt;&lt; <span class="string">"min "</span> &lt;&lt; thermo1.<span class="type">T</span><span class="literal">()</span>.name<span class="literal">()</span></span><br><span class="line">       &lt;&lt; <span class="string">" "</span> &lt;&lt; min(thermo1.<span class="type">T</span><span class="literal">()</span>).<span class="keyword">value</span><span class="literal">()</span> &lt;&lt; endl;</span><br><span class="line"></span><br><span class="line">   thermo2.correct<span class="literal">()</span>;</span><br><span class="line">   <span class="type">Info</span>&lt;&lt; <span class="string">"min "</span> &lt;&lt; thermo2.<span class="type">T</span><span class="literal">()</span>.name<span class="literal">()</span></span><br><span class="line">       &lt;&lt; <span class="string">" "</span> &lt;&lt; min(thermo2.<span class="type">T</span><span class="literal">()</span>).<span class="keyword">value</span><span class="literal">()</span> &lt;&lt; endl;</span><br></pre></td></tr></table></figure></p>
<p>这里， <code>solve</code> 函数值得细说，重点是看 <code>correct()</code> 函数，以及 <code>T()</code> 函数。</p>
<p><code>corretc()</code> 函数指的是定义在 <code>heRhoThermo</code> 类中的 <code>correct()</code> 函数：<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">BasicPsiThermo</span>, class <span class="type">MixtureType</span>&gt;</span><br><span class="line"><span class="type">void</span> <span class="type">Foam</span>::heRhoThermo&lt;<span class="type">BasicPsiThermo</span>, <span class="type">MixtureType</span>&gt;::correct()</span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">if</span> (debug)</span><br><span class="line">    &#123;</span><br><span class="line">        <span class="type">Info</span>&lt;&lt; <span class="string">"entering heRhoThermo&lt;MixtureType&gt;::correct()"</span> &lt;&lt; endl;</span><br><span class="line">    &#125;</span><br><span class="line"></span><br><span class="line">    calculate();</span><br><span class="line"></span><br><span class="line">    <span class="keyword">if</span> (debug)</span><br><span class="line">    &#123;</span><br><span class="line">        <span class="type">Info</span>&lt;&lt; <span class="string">"exiting heRhoThermo&lt;MixtureType&gt;::correct()"</span> &lt;&lt; endl;</span><br><span class="line">    &#125;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>可见， <code>correct()</code> 函数，其实就是对 <code>calculate</code> 函数进行了一次调用而已。<br>看来最核心最关键的就在 <code>calculate</code> 函数中了。在仔细看这个函数之前，先把 <code>T()</code> 的定义看完。 <code>T()</code> 定义在 <code>basicThermo</code> 类中，其作用仅是返回同样定义在 <code>basicThermo</code> 类中定义的数据成员 <code>T_</code> 而已。</p>
<p>下面深入分析一下 <code>heRhoThermo</code> 类中的 <code>calculate</code> 函数，这里再将它列出来一次：<br><figure class="highlight prolog"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br><span class="line">26</span><br><span class="line">27</span><br><span class="line">28</span><br><span class="line">29</span><br><span class="line">30</span><br><span class="line">31</span><br><span class="line">32</span><br><span class="line">33</span><br><span class="line">34</span><br><span class="line">35</span><br><span class="line">36</span><br><span class="line">37</span><br><span class="line">38</span><br><span class="line">39</span><br><span class="line">40</span><br><span class="line">41</span><br><span class="line">42</span><br><span class="line">43</span><br><span class="line">44</span><br><span class="line">45</span><br><span class="line">46</span><br><span class="line">47</span><br><span class="line">48</span><br><span class="line">49</span><br><span class="line">50</span><br><span class="line">51</span><br><span class="line">52</span><br><span class="line">53</span><br><span class="line">54</span><br><span class="line">55</span><br><span class="line">56</span><br><span class="line">57</span><br><span class="line">58</span><br><span class="line">59</span><br><span class="line">60</span><br><span class="line">61</span><br><span class="line">62</span><br><span class="line">63</span><br><span class="line">64</span><br><span class="line">65</span><br><span class="line">66</span><br><span class="line">67</span><br><span class="line">68</span><br><span class="line">69</span><br><span class="line">70</span><br><span class="line">71</span><br><span class="line">72</span><br><span class="line">73</span><br><span class="line">74</span><br><span class="line">75</span><br></pre></td><td class="code"><pre><span class="line"><span class="atom">template</span>&lt;<span class="atom">class</span> <span class="name">BasicPsiThermo</span>, <span class="atom">class</span> <span class="name">MixtureType</span>&gt;</span><br><span class="line"><span class="atom">void</span> <span class="name">Foam</span>::<span class="atom">heRhoThermo</span>&lt;<span class="name">BasicPsiThermo</span>, <span class="name">MixtureType</span>&gt;::<span class="atom">calculate</span>()</span><br><span class="line">&#123;</span><br><span class="line">    <span class="atom">const</span> <span class="atom">scalarField</span>&amp; <span class="atom">hCells</span> = <span class="atom">this</span>-&gt;<span class="atom">he</span>().<span class="atom">internalField</span>();</span><br><span class="line">    <span class="atom">const</span> <span class="atom">scalarField</span>&amp; <span class="atom">pCells</span> = <span class="atom">this</span>-&gt;<span class="atom">p_</span>.<span class="atom">internalField</span>();</span><br><span class="line"></span><br><span class="line">    <span class="atom">scalarField</span>&amp; <span class="name">TCells</span> = <span class="atom">this</span>-&gt;<span class="name">T_</span>.<span class="atom">internalField</span>();</span><br><span class="line">    <span class="atom">scalarField</span>&amp; <span class="atom">psiCells</span> = <span class="atom">this</span>-&gt;<span class="atom">psi_</span>.<span class="atom">internalField</span>();</span><br><span class="line">    <span class="atom">scalarField</span>&amp; <span class="atom">rhoCells</span> = <span class="atom">this</span>-&gt;<span class="atom">rho_</span>.<span class="atom">internalField</span>();</span><br><span class="line">    <span class="atom">scalarField</span>&amp; <span class="atom">muCells</span> = <span class="atom">this</span>-&gt;<span class="atom">mu_</span>.<span class="atom">internalField</span>();</span><br><span class="line">    <span class="atom">scalarField</span>&amp; <span class="atom">alphaCells</span> = <span class="atom">this</span>-&gt;<span class="atom">alpha_</span>.<span class="atom">internalField</span>();</span><br><span class="line"></span><br><span class="line">    <span class="atom">forAll</span>(<span class="name">TCells</span>, <span class="atom">celli</span>)</span><br><span class="line">    &#123;</span><br><span class="line">        <span class="atom">const</span> <span class="atom">typename</span> <span class="name">MixtureType</span>::<span class="atom">thermoType</span>&amp; <span class="atom">mixture_</span> =</span><br><span class="line">            <span class="atom">this</span>-&gt;<span class="atom">cellMixture</span>(<span class="atom">celli</span>);</span><br><span class="line"></span><br><span class="line">        <span class="name">TCells</span>[<span class="atom">celli</span>] = <span class="atom">mixture_</span>.<span class="name">THE</span></span><br><span class="line">        (</span><br><span class="line">            <span class="atom">hCells</span>[<span class="atom">celli</span>],</span><br><span class="line">            <span class="atom">pCells</span>[<span class="atom">celli</span>],</span><br><span class="line">            <span class="name">TCells</span>[<span class="atom">celli</span>]</span><br><span class="line">        );</span><br><span class="line"></span><br><span class="line">        <span class="atom">psiCells</span>[<span class="atom">celli</span>] = <span class="atom">mixture_</span>.<span class="atom">psi</span>(<span class="atom">pCells</span>[<span class="atom">celli</span>], <span class="name">TCells</span>[<span class="atom">celli</span>]);</span><br><span class="line">        <span class="atom">rhoCells</span>[<span class="atom">celli</span>] = <span class="atom">mixture_</span>.<span class="atom">rho</span>(<span class="atom">pCells</span>[<span class="atom">celli</span>], <span class="name">TCells</span>[<span class="atom">celli</span>]);</span><br><span class="line"></span><br><span class="line">        <span class="atom">muCells</span>[<span class="atom">celli</span>] = <span class="atom">mixture_</span>.<span class="atom">mu</span>(<span class="atom">pCells</span>[<span class="atom">celli</span>], <span class="name">TCells</span>[<span class="atom">celli</span>]);</span><br><span class="line">        <span class="atom">alphaCells</span>[<span class="atom">celli</span>] = <span class="atom">mixture_</span>.<span class="atom">alphah</span>(<span class="atom">pCells</span>[<span class="atom">celli</span>], <span class="name">TCells</span>[<span class="atom">celli</span>]);</span><br><span class="line">    &#125;</span><br><span class="line"></span><br><span class="line">    <span class="atom">forAll</span>(<span class="atom">this</span>-&gt;<span class="name">T_</span>.<span class="atom">boundaryField</span>(), <span class="atom">patchi</span>)</span><br><span class="line">    &#123;</span><br><span class="line">        <span class="atom">fvPatchScalarField</span>&amp; <span class="atom">pp</span> = <span class="atom">this</span>-&gt;<span class="atom">p_</span>.<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line">        <span class="atom">fvPatchScalarField</span>&amp; <span class="atom">pT</span> = <span class="atom">this</span>-&gt;<span class="name">T_</span>.<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line">        <span class="atom">fvPatchScalarField</span>&amp; <span class="atom">ppsi</span> = <span class="atom">this</span>-&gt;<span class="atom">psi_</span>.<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line">        <span class="atom">fvPatchScalarField</span>&amp; <span class="atom">prho</span> = <span class="atom">this</span>-&gt;<span class="atom">rho_</span>.<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line"></span><br><span class="line">        <span class="atom">fvPatchScalarField</span>&amp; <span class="atom">ph</span> = <span class="atom">this</span>-&gt;<span class="atom">he</span>().<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line"></span><br><span class="line">        <span class="atom">fvPatchScalarField</span>&amp; <span class="atom">pmu</span> = <span class="atom">this</span>-&gt;<span class="atom">mu_</span>.<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line">        <span class="atom">fvPatchScalarField</span>&amp; <span class="atom">palpha</span> = <span class="atom">this</span>-&gt;<span class="atom">alpha_</span>.<span class="atom">boundaryField</span>()[<span class="atom">patchi</span>];</span><br><span class="line"></span><br><span class="line">        <span class="atom">if</span> (<span class="atom">pT</span>.<span class="atom">fixesValue</span>())</span><br><span class="line">        &#123;</span><br><span class="line">            <span class="atom">forAll</span>(<span class="atom">pT</span>, <span class="atom">facei</span>)</span><br><span class="line">            &#123;</span><br><span class="line">                <span class="atom">const</span> <span class="atom">typename</span> <span class="name">MixtureType</span>::<span class="atom">thermoType</span>&amp; <span class="atom">mixture_</span> =</span><br><span class="line">                    <span class="atom">this</span>-&gt;<span class="atom">patchFaceMixture</span>(<span class="atom">patchi</span>, <span class="atom">facei</span>);</span><br><span class="line"></span><br><span class="line">                <span class="atom">ph</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="name">HE</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line"></span><br><span class="line">                <span class="atom">ppsi</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">psi</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">                <span class="atom">prho</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">rho</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">                <span class="atom">pmu</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">mu</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">                <span class="atom">palpha</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">alphah</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">            &#125;</span><br><span class="line">        &#125;</span><br><span class="line">        <span class="atom">else</span></span><br><span class="line">        &#123;</span><br><span class="line">            <span class="atom">forAll</span>(<span class="atom">pT</span>, <span class="atom">facei</span>)</span><br><span class="line">            &#123;</span><br><span class="line">                <span class="atom">const</span> <span class="atom">typename</span> <span class="name">MixtureType</span>::<span class="atom">thermoType</span>&amp; <span class="atom">mixture_</span> =</span><br><span class="line">                    <span class="atom">this</span>-&gt;<span class="atom">patchFaceMixture</span>(<span class="atom">patchi</span>, <span class="atom">facei</span>);</span><br><span class="line"></span><br><span class="line">                <span class="atom">pT</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="name">THE</span>(<span class="atom">ph</span>[<span class="atom">facei</span>], <span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line"></span><br><span class="line">                <span class="atom">ppsi</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">psi</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">                <span class="atom">prho</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">rho</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">                <span class="atom">pmu</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">mu</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">                <span class="atom">palpha</span>[<span class="atom">facei</span>] = <span class="atom">mixture_</span>.<span class="atom">alphah</span>(<span class="atom">pp</span>[<span class="atom">facei</span>], <span class="atom">pT</span>[<span class="atom">facei</span>]);</span><br><span class="line">            &#125;</span><br><span class="line">        &#125;</span><br><span class="line">    &#125;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这个函数是在对几个热物理相关的量来进行更新，先更新内部场，再更新边界值。一个一个来看：</p>
<ul>
<li><p>he<br>he 前面讲了，这里需要注意的是其边界值的更新。由于 <code>he</code> 没有IO，其内部场量通过求解能量方程来更新，边界则需要根据情况特殊处理。两种情况，一种是设定了边界的温度值(pT.fixesValue())，这时需要更新边界上的 <code>he</code> 值： <code>ph[facei] = mixture_.HE(pp[facei], pT[facei]);</code> 否则，则边界上的 <code>he</code> 不需要特殊地更新。</p>
</li>
<li><p>psi<br>这个直接调用两参数的 <code>psi</code> 函数来更新，这个函数的定义在状态方程里，以 <code>perfaceGas</code> 为例，</p>
<figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;<span class="keyword">class</span> Specie&gt;</span><br><span class="line"><span class="keyword">inline</span> Foam::scalar Foam::perfectGas&lt;Specie&gt;::psi(scalar, scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> <span class="number">1.0</span>/(<span class="keyword">this</span>-&gt;R()*T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure>
</li>
</ul>
<p><code>psi</code> 是压缩因子，返回 $\frac{1}{RT}$。</p>
<ul>
<li>rho<br>这个调用的是两参数的 <code>rho</code> 函数，定义在状态方程类中，用于密度的更新，同样以 <code>perfaceGas</code> 为例，<figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Specie</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar <span class="type">Foam</span>::perfectGas&lt;<span class="type">Specie</span>&gt;::rho(scalar p, scalar T) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> p/(this-&gt;R()*T);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure>
</li>
</ul>
<p>返回 $\frac{p}{RT}$。</p>
<ul>
<li>mu<br>这个调用的是两参数的 <code>mu</code> 函数，其定义在 transport 类中，以 <code>constTransport</code> 为例，这个返回的是场量的层流粘度<figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar <span class="type">Foam</span>::constTransport&lt;<span class="type">Thermo</span>&gt;::mu</span><br><span class="line">(</span><br><span class="line">     <span class="keyword">const</span> scalar p,</span><br><span class="line">     <span class="keyword">const</span> scalar T</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">     <span class="keyword">return</span> mu_; // 常量</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure>
</li>
</ul>
<p>alphah 上面说过了，不再重复。</p>
<p>最后，最复杂的就是温度的更新了</p>
<ul>
<li>T<br>温度的更新，调用的是三参数的 <code>THE</code> 函数，这个函数定义在 <code>species::thermo</code> 类中，<figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar <span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::<span class="type">THE</span></span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> scalar he,</span><br><span class="line">    <span class="keyword">const</span> scalar p,</span><br><span class="line">    <span class="keyword">const</span> scalar <span class="type">T0</span></span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> <span class="type">Type</span>&lt;thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt; &gt;::<span class="type">THE</span>(*this, he, p, <span class="type">T0</span>);</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure>
</li>
</ul>
<p>这里，调用的是 <code>energy variable</code> 类的 <code>THE</code> 函数，以 <code>sensibleInternalEnergy</code> 为例，<br><figure class="highlight stata"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">scalar</span> THE</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> Thermo&amp; thermo,</span><br><span class="line">    <span class="keyword">const</span> <span class="keyword">scalar</span> <span class="keyword">e</span>,</span><br><span class="line">    <span class="keyword">const</span> <span class="keyword">scalar</span> p,</span><br><span class="line">    <span class="keyword">const</span> <span class="keyword">scalar</span> T0</span><br><span class="line"> ) <span class="keyword">const</span></span><br><span class="line"> &#123;</span><br><span class="line">    <span class="keyword">return</span> thermo.<span class="keyword">TEs</span>(<span class="keyword">e</span>, p, T0);</span><br><span class="line"> &#125;</span><br></pre></td></tr></table></figure></p>
<p>可见，对于 <code>sensibleInternalEnergy</code> ， <code>THE</code> 函数实际上返回的是 <code>species::thermo</code> 类的 <code>TEs</code> 函数。<br><figure class="highlight nimrod"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;class <span class="type">Thermo</span>, <span class="keyword">template</span>&lt;class&gt; class <span class="type">Type</span>&gt;</span><br><span class="line">inline <span class="type">Foam</span>::scalar <span class="type">Foam</span>::species::thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::<span class="type">TEs</span></span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> scalar es,</span><br><span class="line">    <span class="keyword">const</span> scalar p,</span><br><span class="line">    <span class="keyword">const</span> scalar <span class="type">T0</span></span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> T</span><br><span class="line">    (</span><br><span class="line">        es,</span><br><span class="line">        p,</span><br><span class="line">        <span class="type">T0</span>,</span><br><span class="line">        &amp;thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::<span class="type">Es</span>,</span><br><span class="line">        &amp;thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::<span class="type">Cv</span>,</span><br><span class="line">        &amp;thermo&lt;<span class="type">Thermo</span>, <span class="type">Type</span>&gt;::limit</span><br><span class="line">    );</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这里，终于来到了这个六参数的 <code>T</code> 函数：<br><figure class="highlight stata"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br><span class="line">26</span><br><span class="line">27</span><br><span class="line">28</span><br><span class="line">29</span><br><span class="line">30</span><br><span class="line">31</span><br><span class="line">32</span><br><span class="line">33</span><br><span class="line">34</span><br><span class="line">35</span><br><span class="line">36</span><br><span class="line">37</span><br><span class="line">38</span><br><span class="line">39</span><br><span class="line">40</span><br><span class="line">41</span><br><span class="line">42</span><br><span class="line">43</span><br><span class="line">44</span><br><span class="line">45</span><br><span class="line">46</span><br><span class="line">47</span><br><span class="line">48</span><br><span class="line">49</span><br><span class="line">50</span><br><span class="line">51</span><br><span class="line">52</span><br><span class="line">53</span><br><span class="line">54</span><br><span class="line">55</span><br><span class="line">56</span><br></pre></td><td class="code"><pre><span class="line"><span class="comment">// 声明</span></span><br><span class="line"> inline <span class="keyword">scalar</span> T</span><br><span class="line">        (</span><br><span class="line">            <span class="keyword">scalar</span> f,</span><br><span class="line">            <span class="keyword">scalar</span> p,</span><br><span class="line">            <span class="keyword">scalar</span> T0,</span><br><span class="line">            <span class="keyword">scalar</span> (thermo::*F)(<span class="keyword">const</span> <span class="keyword">scalar</span>, <span class="keyword">const</span> <span class="keyword">scalar</span>) <span class="keyword">const</span>,</span><br><span class="line">            <span class="keyword">scalar</span> (thermo::*dFdT)(<span class="keyword">const</span> <span class="keyword">scalar</span>, <span class="keyword">const</span> <span class="keyword">scalar</span>) <span class="keyword">const</span>,</span><br><span class="line">            <span class="keyword">scalar</span> (thermo::*limit)(<span class="keyword">const</span> <span class="keyword">scalar</span>) <span class="keyword">const</span></span><br><span class="line">        ) <span class="keyword">const</span>;</span><br><span class="line"></span><br><span class="line"><span class="comment">//实现</span></span><br><span class="line">template&lt;<span class="keyword">class</span> Thermo, template&lt;<span class="keyword">class</span>&gt; <span class="keyword">class</span> <span class="keyword">Type</span>&gt;</span><br><span class="line">inline Foam::<span class="keyword">scalar</span> Foam::species::thermo&lt;Thermo, <span class="keyword">Type</span>&gt;::T</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">scalar</span> f,</span><br><span class="line">    <span class="keyword">scalar</span> p,</span><br><span class="line">    <span class="keyword">scalar</span> T0,</span><br><span class="line">    <span class="keyword">scalar</span> (thermo&lt;Thermo, <span class="keyword">Type</span>&gt;::*F)(<span class="keyword">const</span> <span class="keyword">scalar</span>, <span class="keyword">const</span> <span class="keyword">scalar</span>) <span class="keyword">const</span>,</span><br><span class="line">    <span class="keyword">scalar</span> (thermo&lt;Thermo, <span class="keyword">Type</span>&gt;::*dFdT)(<span class="keyword">const</span> <span class="keyword">scalar</span>, <span class="keyword">const</span> <span class="keyword">scalar</span>)</span><br><span class="line">        <span class="keyword">const</span>,</span><br><span class="line">    <span class="keyword">scalar</span> (thermo&lt;Thermo, <span class="keyword">Type</span>&gt;::*limit)(<span class="keyword">const</span> <span class="keyword">scalar</span>) <span class="keyword">const</span></span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">scalar</span> <span class="keyword">Test</span> = T0;</span><br><span class="line">    <span class="keyword">scalar</span> Tnew = T0;</span><br><span class="line">    <span class="keyword">scalar</span> Ttol = T0*tol_;</span><br><span class="line">    int    iter = 0;</span><br><span class="line"></span><br><span class="line">    <span class="keyword">do</span></span><br><span class="line">    &#123;</span><br><span class="line">        <span class="keyword">Test</span> = Tnew;</span><br><span class="line">        Tnew =</span><br><span class="line">            (this-&gt;*limit)</span><br><span class="line">            (<span class="keyword">Test</span> - ((this-&gt;*F)(p, <span class="keyword">Test</span>) - f)/(this-&gt;*dFdT)(p, <span class="keyword">Test</span>));</span><br><span class="line"></span><br><span class="line">        <span class="keyword">if</span> (iter++ &gt; maxIter_)</span><br><span class="line">        &#123;</span><br><span class="line">            FatalErrorIn</span><br><span class="line">            (</span><br><span class="line">                <span class="string">"thermo&lt;Thermo, Type&gt;::T(scalar f, scalar T0, "</span></span><br><span class="line">                <span class="string">"scalar (thermo&lt;Thermo, Type&gt;::*F)"</span></span><br><span class="line">                <span class="string">"(const scalar) const, "</span></span><br><span class="line">                <span class="string">"scalar (thermo&lt;Thermo, Type&gt;::*dFdT)"</span></span><br><span class="line">                <span class="string">"(const scalar) const, "</span></span><br><span class="line">                <span class="string">"scalar (thermo&lt;Thermo, Type&gt;::*limit)"</span></span><br><span class="line">                <span class="string">"(const scalar) const"</span></span><br><span class="line">                <span class="string">") const"</span></span><br><span class="line">            )   &lt;&lt; <span class="string">"Maximum number of iterations exceeded"</span></span><br><span class="line">                &lt;&lt; abort(FatalError);</span><br><span class="line">        &#125;</span><br><span class="line"></span><br><span class="line">    &#125; <span class="keyword">while</span> (mag(Tnew - <span class="keyword">Test</span>) &gt; Ttol);</span><br><span class="line"></span><br><span class="line">    <span class="keyword">return</span> Tnew;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>这个函数，前三个参数是普通的 scalar 类型变量，后三个参数，是函数指针，并且都指向当前类 <code>species::thermo</code> 的成员函数。以 <code>TEs</code> 为例，后三个参数分别代入的是 <code>Es</code> ， <code>Cv</code> 以及 <code>limit</code> 三个函数。 <code>Es</code> 和 <code>Cv</code> 前面都看过了， <code>limit</code> 定义在 thermo 类中，以 <code>hConst</code> 为例，<br><figure class="highlight cpp"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">template</span>&lt;<span class="keyword">class</span> EquationOfState&gt;</span><br><span class="line"><span class="keyword">inline</span> Foam::scalar Foam::hConstThermo&lt;EquationOfState&gt;::limit</span><br><span class="line">(</span><br><span class="line">    <span class="keyword">const</span> scalar T</span><br><span class="line">) <span class="keyword">const</span></span><br><span class="line">&#123;</span><br><span class="line">    <span class="keyword">return</span> T;</span><br><span class="line">&#125;</span><br></pre></td></tr></table></figure></p>
<p>直接返回温度 <code>T</code> 。事实上，除了 <code>janaf</code> 模型，其他的都是返回 <code>T</code> 。 <code>janaf</code> 模型中， 如果温度没有超出 [Tlow,Thigh]，则会出来警告信息，并且，若 <code>T &lt; Tlow</code> 则返回 <code>Tlow</code>，而 <code>T &gt; Thigh</code> 时，则返回 <code>Thigh</code>。</p>
<p>下面仔细来分析六参数 <code>T</code> 函数的核心部分。经过摸索，发现这个其实是一个牛顿迭代的过程，目的是根据 <code>Es</code> 函数，从内能 <code>es</code> 来计算温度 <code>T</code>，即求解 $E_s(p,T) - E_s = 0$ 。令 $F(T)= E_s(p,T) - E_s $，则牛顿迭代法的递推公式为<br>$$<br>T_{New} = T_{old} - \dfrac{F(T_{old})}{F\prime(T_{old})} = T_{old} - \dfrac{E_s(p, T_{old)} - E_s}{dE_s(p,T)/dT |_{T=T_{old}}}<br>$$<br>对于 <code>sensibleInternalEnergy</code> ，$dE_s(p,T)/dT = C_v(p,T)$<br>所以最终得到递推公式为<br>$$<br>T_{New} = T_{old} - \dfrac{E_s(p, T_{old)} - E_s}{C_v(p, T_{old})}<br>$$<br>这里设置了最大迭代次数为 100，超过将报那个涉及到能量的模拟中最容易见到的崩溃信息：”Maximum number of iterations exceeded” 。</p>
<p>当能量变量是焓时，则 $E_s$ 要换成 $H_s$， $C_v$ 要换成 $C_p$ 。 </p>
<p>至此便分析完了一个具体的能量方程实例。</p>

      
    </div>
    <footer class="article-footer">
      <a data-url="http://xiaopingqiu.github.io/2016/06/25/thermophysics4/" data-id="cj2lay31h001hkkmbzazmdbk9" class="article-share-link">Share</a>
      
      
  <ul class="article-tag-list"><li class="article-tag-list-item"><a class="article-tag-list-link" href="/tags/Code-Explained/">Code Explained</a></li><li class="article-tag-list-item"><a class="article-tag-list-link" href="/tags/thermophysicalModels/">thermophysicalModels</a></li></ul>

    </footer>
  </div>
  
    
<nav id="article-nav">
  
    <a href="/2016/06/25/thermophysics5/" id="article-nav-newer" class="article-nav-link-wrap">
      <strong class="article-nav-caption">Newer</strong>
      <div class="article-nav-title">
        
          OpenFOAM 中的热物理类之添加新模型
        
      </div>
    </a>
  
  
    <a href="/2016/06/25/thermophysics3/" id="article-nav-older" class="article-nav-link-wrap">
      <strong class="article-nav-caption">Older</strong>
      <div class="article-nav-title">OpenFOAM 中的热物理类之继承派生关系</div>
    </a>
  
</nav>

  
</article>


<section id="comments">
  <!-- 多说评论框 start -->
  <div class="ds-thread" data-thread-key="post-thermophysics4" data-title="一个具体能量方程的解析" data-url="http://xiaopingqiu.github.io/2016/06/25/thermophysics4/"></div>
  <!-- 多说评论框 end -->
  <!-- 多说公共JS代码 start (一个网页只需插入一次) -->
  <script type="text/javascript">
  var duoshuoQuery = {short_name:"xiaopingqiu"};
	(function() {
		var ds = document.createElement('script');
		ds.type = 'text/javascript';ds.async = true;
		ds.src = (document.location.protocol == 'https:' ? 'https:' : 'http:') + '//static.duoshuo.com/embed.js';
		ds.charset = 'UTF-8';
		(document.getElementsByTagName('head')[0] 
		 || document.getElementsByTagName('body')[0]).appendChild(ds);
	})();
	</script>
  <!-- 多说公共JS代码 end -->
</section>

</section>
        
          <aside id="sidebar">
  
    
  <div class="widget-wrap">
    <h3 class="widget-title">Categories</h3>
    <div class="widget">
      <ul class="category-list"><li class="category-list-item"><a class="category-list-link" href="/categories/C/">C++</a><span class="category-list-count">2</span></li><li class="category-list-item"><a class="category-list-link" href="/categories/DEM/">DEM</a><span class="category-list-count">1</span></li><li class="category-list-item"><a class="category-list-link" href="/categories/OpenFOAM/">OpenFOAM</a><span class="category-list-count">44</span></li><li class="category-list-item"><a class="category-list-link" href="/categories/Paraview/">Paraview</a><span class="category-list-count">5</span></li><li class="category-list-item"><a class="category-list-link" href="/categories/swak4Foam/">swak4Foam</a><span class="category-list-count">1</span></li><li class="category-list-item"><a class="category-list-link" href="/categories/test/">test</a><span class="category-list-count">2</span></li><li class="category-list-item"><a class="category-list-link" href="/categories/vim/">vim</a><span class="category-list-count">1</span></li></ul>
    </div>
  </div>

  
    
  <div class="widget-wrap">
    <h3 class="widget-title">Tags</h3>
    <div class="widget">
      <ul class="tag-list"><li class="tag-list-item"><a class="tag-list-link" href="/tags/Boundary-conditions/">Boundary conditions</a><span class="tag-list-count">6</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/C/">C++</a><span class="tag-list-count">2</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/CentOS/">CentOS</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/Code-Explained/">Code Explained</a><span class="tag-list-count">29</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/LES/">LES</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/LIGGGHTS/">LIGGGHTS</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/ODE/">ODE</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/OpenFOAM/">OpenFOAM</a><span class="tag-list-count">20</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/Postprocessing/">Postprocessing</a><span class="tag-list-count">9</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/Preprocessing/">Preprocessing</a><span class="tag-list-count">2</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/RTS/">RTS</a><span class="tag-list-count">3</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/TIL/">TIL</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/Windows/">Windows</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/fvOptions/">fvOptions</a><span class="tag-list-count">2</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/groovyBC/">groovyBC</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/paraview/">paraview</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/test/">test</a><span class="tag-list-count">2</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/thermophysicalModels/">thermophysicalModels</a><span class="tag-list-count">5</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/turbulence-model/">turbulence model</a><span class="tag-list-count">7</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/vim/">vim</a><span class="tag-list-count">1</span></li><li class="tag-list-item"><a class="tag-list-link" href="/tags/wall-functions/">wall functions</a><span class="tag-list-count">4</span></li></ul>
    </div>
  </div>

  
    
  <div class="widget-wrap">
    <h3 class="widget-title">Tag Cloud</h3>
    <div class="widget tagcloud">
      <a href="/tags/Boundary-conditions/" style="font-size: 15.56px;">Boundary conditions</a><a href="/tags/C/" style="font-size: 11.11px;">C++</a><a href="/tags/CentOS/" style="font-size: 10px;">CentOS</a><a href="/tags/Code-Explained/" style="font-size: 20px;">Code Explained</a><a href="/tags/LES/" style="font-size: 10px;">LES</a><a href="/tags/LIGGGHTS/" style="font-size: 10px;">LIGGGHTS</a><a href="/tags/ODE/" style="font-size: 10px;">ODE</a><a href="/tags/OpenFOAM/" style="font-size: 18.89px;">OpenFOAM</a><a href="/tags/Postprocessing/" style="font-size: 17.78px;">Postprocessing</a><a href="/tags/Preprocessing/" style="font-size: 11.11px;">Preprocessing</a><a href="/tags/RTS/" style="font-size: 12.22px;">RTS</a><a href="/tags/TIL/" style="font-size: 10px;">TIL</a><a href="/tags/Windows/" style="font-size: 10px;">Windows</a><a href="/tags/fvOptions/" style="font-size: 11.11px;">fvOptions</a><a href="/tags/groovyBC/" style="font-size: 10px;">groovyBC</a><a href="/tags/paraview/" style="font-size: 10px;">paraview</a><a href="/tags/test/" style="font-size: 11.11px;">test</a><a href="/tags/thermophysicalModels/" style="font-size: 14.44px;">thermophysicalModels</a><a href="/tags/turbulence-model/" style="font-size: 16.67px;">turbulence model</a><a href="/tags/vim/" style="font-size: 10px;">vim</a><a href="/tags/wall-functions/" style="font-size: 13.33px;">wall functions</a>
    </div>
  </div>

  
    
  <div class="widget-wrap">
    <h3 class="widget-title">Archives</h3>
    <div class="widget">
      <ul class="archive-list"><li class="archive-list-item"><a class="archive-list-link" href="/archives/2017/05/">五月 2017</a><span class="archive-list-count">1</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2016/08/">八月 2016</a><span class="archive-list-count">8</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2016/06/">六月 2016</a><span class="archive-list-count">5</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2016/05/">五月 2016</a><span class="archive-list-count">3</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2016/04/">四月 2016</a><span class="archive-list-count">12</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2016/03/">三月 2016</a><span class="archive-list-count">6</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/12/">十二月 2015</a><span class="archive-list-count">1</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/11/">十一月 2015</a><span class="archive-list-count">2</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/10/">十月 2015</a><span class="archive-list-count">1</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/09/">九月 2015</a><span class="archive-list-count">6</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/08/">八月 2015</a><span class="archive-list-count">1</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/06/">六月 2015</a><span class="archive-list-count">2</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/05/">五月 2015</a><span class="archive-list-count">6</span></li><li class="archive-list-item"><a class="archive-list-link" href="/archives/2015/04/">四月 2015</a><span class="archive-list-count">2</span></li></ul>
    </div>
  </div>

  
    
  <div class="widget-wrap">
    <h3 class="widget-title">Recents</h3>
    <div class="widget">
      <ul>
        
          <li>
            <a href="/2017/05/12/duoshuoAnnouncement/">多说评论系统将停止提供服务</a>
          </li>
        
          <li>
            <a href="/2016/08/27/ParaviewScritps/">Paraview 脚本一例</a>
          </li>
        
          <li>
            <a href="/2016/08/27/ParaviewCamera/">Paraview 中有关 Camera 的操作两例</a>
          </li>
        
          <li>
            <a href="/2016/08/27/ParaviewCustomFilter/">Paraview 中创建 Custom Filter</a>
          </li>
        
          <li>
            <a href="/2016/08/27/ParaviewStreamLineOnSlice/">在 Paraview 中画截面上的流线</a>
          </li>
        
      </ul>
    </div>
  </div>

  
</aside>
        
      </div>
      <footer id="footer">
  
  <div class="outer">
    <div id="footer-info" class="inner">
      &copy; 2017 Giskard Q.<br>
      Powered by <a href="http://hexo.io/" target="_blank">Hexo</a>
    </div>
  </div>
</footer>

<script src="//dn-lbstatics.qbox.me/lbservice/busuanzi/2.0/busuanzi.mini.js"/></script>

    </div>
    <nav id="mobile-nav">
  
    <a href="/" class="mobile-nav-link">Home</a>
  
    <a href="/archives" class="mobile-nav-link">Archives</a>
  
    <a href="/atom.xml" class="mobile-nav-link">Rss</a>
  
    <a href="/about" class="mobile-nav-link">About</a>
  
</nav>
    

<script src="//ajax.googleapis.com/ajax/libs/jquery/2.0.3/jquery.min.js"></script>


  <link rel="stylesheet" href="/fancybox/jquery.fancybox.css" type="text/css">
  <script src="/fancybox/jquery.fancybox.pack.js" type="text/javascript"></script>


<script src="/js/script.js" type="text/javascript"></script>

<script>
  (function(i,s,o,g,r,a,m){i['GoogleAnalyticsObject']=r;i[r]=i[r]||function(){
  (i[r].q=i[r].q||[]).push(arguments)},i[r].l=1*new Date();a=s.createElement(o),
  m=s.getElementsByTagName(o)[0];a.async=1;a.src=g;m.parentNode.insertBefore(a,m)
  })(window,document,'script','//www.google-analytics.com/analytics.js','ga');

  ga('create', 'UA-62501686-1', 'auto');
  ga('send', 'pageview');

</script>

  </div>
<!-- mathjax config similar to math.stackexchange -->

<script type="text/x-mathjax-config">
  MathJax.Hub.Config({
    tex2jax: {
      inlineMath: [ ['$','$'], ["\\(","\\)"] ],
      processEscapes: true
    }
  });
</script>

<script type="text/x-mathjax-config">
    MathJax.Hub.Config({
      tex2jax: {
        skipTags: ['script', 'noscript', 'style', 'textarea', 'pre', 'code']
      }
    });
</script>

<script type="text/x-mathjax-config">
    MathJax.Hub.Queue(function() {
        var all = MathJax.Hub.getAllJax(), i;
        for(i=0; i < all.length; i += 1) {
            all[i].SourceElement().parentNode.className += ' has-jax';
        }
    });
</script>

<script type="text/javascript" src="http://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML">
</script>
</body>
</html>